#!/bin/bash -e

function info() {
echo Usage: `basename $0` 'in.vcf'
exit 1
}

while getopts  ":p:f:" opt; do
	case  $opt  in
		p) out_prefix=$OPTARG;;
		f) suffix=$OPTARG;;
		*) info;;
	esac
done
shift $(($OPTIND - 1))


if [ $# -lt 1 ]; then info; fi

. /mnt/ilustre/app/medical/tools/script/var.sh

# $java_run/snpsift vcfCheck $1

# exp="(1000g2015aug_all < 0.005) & (!(Func.refGene = 'intronic')) & (!(ExonicFunc.refGene = 'synonymous_SNV')) & isHet(GEN[0]) & isHet(GEN[1]) & isHom(GEN[2])"
# exp1="isHet(GEN[0]) & isHet(GEN[1]) & isHom(GEN[2])"
# echo $exp1
# cat $1 | $java_run/snpsift filter "(1000g2015aug_all < 0.01) & (!(Func.refGene = 'intronic')) & (!(ExonicFunc.refGene = 'synonymous_SNV')) & isHet(GEN[0]) & isHet(GEN[1]) & isHom(GEN[2])" > $out_prefix.filt.vcf

cat $1 | $java_run/snpsift filter "(GEN[?].AD[?] > 2)" > $out_prefix.filt.dp.vcf

cat $out_prefix.filt.dp.vcf | $java_run/snpsift filter "(! (1000g2015aug_all > 0.01))" > $out_prefix.filt.af.vcf

cat $out_prefix.filt.af.vcf |$java_run/snpsift filter "(isHet(GEN[0]) & (isHet(GEN[1])) & (isVariant(GEN[2])))" > $out_prefix.filt.gt.vcf

cat $out_prefix.filt.gt.vcf | $java_run/snpsift filter "(! (ExonicFunc.refGene = 'synonymous_SNV'))" > $out_prefix.filt.exon_func.vcf

cat $out_prefix.filt.exon_func.vcf | $java_run/snpsift filter "(! (Func.refGene = 'intronic'))" > $out_prefix.filt.no_intro.vcf


/mnt/ilustre/app/medical/tools/tools/java -jar /mnt/ilustre/app/medical/tools/jars/snpsift extractFields $out_prefix.filt.no_intro.vcf "Gene.refGene" |sed '1d'|sort|uniq -c|sed 's/^  *//; s/ /\t/g'|sort -k1nr > $out_prefix.genes.txt

. $cmd_done